Interplay of orbital and spin ordering in the iron pnictides 
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A number of recent experiments exhibit electronic anisotropy in the iron pnictides, and there 
is a growing body of experimental evidence that its origin is related to orbital ordering in Fe d xz 
and d yz orbitals. We examine this problem in the parent compounds of the iron pnictides by a 
combination of ab initio band theory calculations, phenomenological Ginzburg-Landau theory of 
coupled orbital and magnetic order parameters, and a microscopic mean-field study of the Kugel- 
Khomskii model. We find that orbital ordering is sufficient to explain a number of key experimental 
observations, in particular linear correlation between the orthorhombic lattice distortion and the 
magnetic ordered moment. We predict that orbital polarization should scale as a square of magnetic 
moment close to Tat. Mediated by orbital polarization, the effective spin-spin exchange interactions 
develop anisotropy in the ab-plane, in accord with recent neutron scattering measurements. 
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Proximity of the superconducting phase to antifcr- 
romagnctic (AFM) order in the phase diagram of the 
newly discovered iron-based superconductors, similar to 
the cuprates family, has fueled renewed interest in the 
role of magnetism in unconventional superconductivity. 
Understanding the nature of magnetic order in the iron 
pnictides is an important aspect of these materials, and 
has received much attention [1-6]. Intriguingly, the mag- 
netic order is preceded or coincides with the tetragonal 
to orthorhombic structural transition [1, 2, 4], breaking 
the d lattice symmetry. A phase with spontaneously 
broken rotational symmetry, often referred to as nematic 
phase, has also been recently reported in the pseudogap 
phase of high-temperature cuprate superconductors [7, 8] 
raising the question of its role in superconductivity (SC). 
In the iron pnictides, the nematicity was observed in neu- 
tron scattering [5] and scanning tunneling microscopy [9] 
below the magnetic ordering temperature Tm- It is also 
seen in the resistivity anisotropy measured in detwinned 
Ba(Fei_ x Co a; )2As2 crystals [10, 11] even above Tm- A 
recent observation of highly anisotropic spin-wave dis- 
persion above Tm in the parent compound BaFe2As2 was 
attributed to spontaneous "spin nematicity" [12]. Simi- 
lar anisotropy was also found inside the superconducting 
phase of optimally doped Ba(Fei_ x Co :c )2As2 [13]. 

These observations raise the question of the role of 
broken C4 symmetry and its relation to the magneto- 
structural phase transition in the iron pnictides. Two 
distinct scenarios have been proposed theoretically. The 
first concept, "spin nematic", is based on the sponta- 
neously broken Z2 Ising symmetry between two collinear 
ordering wave-vectors Q a = (n, 0) and Qb = (0, n) cho- 
sen by "order from disorder" mechanism [14-16]. It was 
first proposed by Chandra, Coleman and Larkin that 
such Ising spin symmetry can be spontaneously broken at 
a temperature higher than the magnetic ordering Tm [17], 
where the ordered moment (Mi) = 0, but the Ising vari- 
able a oc (Ma ■ M B ) = ±1 acquires an expectation value 



( A and B are two sublattices). This mechanism, although 
breaks the lattice 64 symmetry, does not however imply a 
broken SU(2) spin-rotational symmetry. In this respect, 
it differs from the term "nematic" originally used in stud- 
ies of liquid crystals and could be perhaps more aptly 
referred to as "spin-Ising ordering", following Rcf. 17. 

The second scenario proposes that the broken or- 
thorhombic symmetry stems from unequal population of 
the d xz and d yz orbitals, resulting in the so-called ferro- 
orbital ordering [18-23]. The order parameter here is the 
orbital polarization P = (n xz — n yz ), which explicitly 
breaks the C4 symmetry. The state with P =/= could 
be called an "electron nematic" and closely corresponds 
to the nematic order proposed for the bilayer ruthenate 
Sr 3 Ru 2 07 [24] under applied magnetic field. The orbital 
ordering has gained support from the recent quadrupolc 
resonance [25] and ARPES measurements [26]. 

Despite very different physical origin, both mecha- 
nisms allow for linear coupling between the structural 
orthorhombic distortion, 5 = and the respective or- 
der parameter, and cannot be distinguished on symmetry 
grounds alone. In this work, we report the combined 
study based on density functional theory, Ginzburg- 
Landau expansion and microscopic theory, all lending 
strong support to the orbital ordering in the pnictides. 

The ab initio density functional theory (DFT) calcu- 
lations have been performed on the parent compound 
BaFe2As2 using the full-potential augmented plane-wave 
basis, as implemented in the WIEN2k code [27], with 
the generalized gradient approximation (GGA) for the 
exchange-correlation functional [28]. Since we are in- 
terested in the spontaneously broken C 4 symmetry, we 
have adopted the experimental tetragonal structure [29], 
to avoid the effects of the orthorhombic lattice distor- 
tion which explicitly breaks the 90° rotation symme- 
try. For the k-point sampling, the 6x6x3 mesh was 
used in the reduced part of the Brillouin zone corre- 
sponding to the 4-atom antiferromagnetic unit cell (4 
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Magnetization in xz and yz orbitals 

FIG. 1: Orbital polarization P — n xz — n yz as function of 
(partial) ordered moment M = | (n X z1 — n X zl + n yz f — n yz i). 
The line is a quadratic fit to the data: P — Po + cpM 2 . The 
inset shows P plotted vs. total magnetization Mtot on the Fe 
site and follows a similar quadratic dependence. 

f.u./cell). The magnetic structure is collinear with or- 
dering wave-vector (ir, 0) or (0, it) in the 1-Fe unit cell 
notation. The calculated ordered moment gM toi =1.91 
fiB is more than twice larger than the measured value 
gM exp — 0.91^s, a known problem of the DFT in 
the iron pnictidcs [30]. One can attempt to include 
the effect of Coulomb repulsion in the LDA+U frame- 
work [31], however usually this exacerbates the prob- 
lem. Recently, it was shown however [32, 33] that a 
low-spin configuration of Fe can be stabilized within 
LDA+J7 approach, in agreement with the tight-binding 
calculations [34]. We indeed find a low- moment solution 
with gMi ow = 0.57/is/Fe within the so-called "around 
mean field" GGA+U scheme [35], using the values of 
U = 2.7 eV and Hund's J = 0.79 eV calculated from a 
constrained- RPA scheme [36] for LaFeAsO [37]. The low- 
moment state is stabilized by a significant energy 0.22 eV 
per Fe atom compared to the high-moment state. 

In both the GGA and GGA+U calculations, we find 
a non-zero value of orbital ordering: Pqga — 0.19 and 
-Pgga+u = 0.30 in the low-moment (tt, 0) phase. The 
GGA calculated polarization agrees well with the value 
reported by the recent ARPES measurement [26] on 
Ba(Fei_ x Coa;)2As2 . In order to elucidate the depen- 
dence of orbital polarization P on the size of the ordered 
moment M, we have performed a series of calculations 
effectively suppressing the spin moment by introducing 
an orbital potential that couples to L z [38]. Shown in 
Fig. 1 is the resulting orbital polarization P = n xz — n yz 
plotted vs. ordered moment in Fe d xz and d yz orbitals 
M = ^(n I2 | — n xz ^ + n yz ^ — n yz ±). The data points can 
be fitted by a quadratic dependence P = Po + cpM 2 , 
which also holds as a function of the total ordered mo- 



ment M to t = \{ n ^ — n i) on Fe site (inset of Fig. 1). 

Phcnomcnologically, the inteplay between orbital and 
spin polarization can be captured by an effective 
Ginzburg-Landau theory with two order parameters: 

F[M, P] = (aP 2 + ^P 4 - wP5^ - 1 P 2 M 2 

+ (rM 2 + |m 4 ) + . . . (1) 

where r = T j^ N * s t ne reduced temperature of the AFM 
transition and S is the orthorhombic lattice distortion. 
We consider the situation below the orbital ordering tran- 
sition temperature, T < To, so that the reduced temper- 
ature a — T yJ° < 0. The choice of the biquadratic cou- 
pling —jP 2 M 2 is warranted by the microscopic Kugel- 
Khomskii theory (see below). Let us first consider the 
case (5 = 0. The saddle-point solution of Eq. (1) yields: 

M= / 2(-r + 7 ps) ifr _ 7p2<Q (2) 
V v 

P = ] j^^J+~^M 2 ^ « P x + c P M 2 + 0{M 4 l?) 

/ \-l/2 

with the coefficient c P {T) = ^P (l + j^yM 2 J and 

Pi (T) = Po ^Jl + j^Mq expressed through the mean- field 

values M (T) = yj-2r/v and P (T) = yj-2a/u in the 
absence of coupling between the two order parameters. 

We see that the Landau theory (3) naturally ac- 
counts for quadratic dependence of orbital polarization 
on the size of the ordered moment, as obtained ear- 
lier from ab initio calculations (Fig. 1). It also fol- 
lows from Eq. (2) that the magnetic ordering tempera- 
ture becomes higher due to coupling to orbital ordering: 
T* N = T N {1 + 2\a*\j/tt) with \a*\ = (T - T* N )/T Q > 0. 

Let us now consider the case of non-zero orthorhom- 
bic distortion S — (a — b)/(a + b). Due to coupling 
of distortion to the orbital ordering in Eq. (1), we ob- 
tain wS = -2aP + uP 3 - 2-fPM 2 . This results in the 
leading linear contribution to orbital order parameter 
P ~ 6 + 0(8 3 ) for temperatures not too far from To- 
It follows from Eq. (2) that the magnetization 

M 2 = — + 2 ^5 2 + 0(S 4 ). (4) 
v v 

In particular near Tjy, where the first term is negligible, 
magnetization scales linearly with orthorhombic distor- 
tion, M (S) oc 5. This prediction of the Landau theory 
is entirely consistent with the neutron scattering mea- 
surements on La(Oi_ x F x )FeAs [1], SrFe2As2 [3, 39] and 
Ca(Fei_ :c Co x ) 2 As 2 [40]. 

Note that instead of biquadratic term — jP 2 M 2 , a cou- 
pling of the form —kPM 2 is also allowed by symmetry, 
however the analysis shows that it would result in the 
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quadratic dependence of structural distortion on mag- 
netic ordered moment, wS = —2aP + uP 3 — kM 2 = 
d + C 2 M 2 + O (M 4 ) . This would clearly contradict the 
aforementioned neutron scattering measurements. 

Intriguingly, the ab initio calculations find that orbital 
polarization is split unequally between the two spin com- 
ponents. Both in GGA and GGA+J7, the largest contri- 
bution comes from the minority spin component. In other 
words, the expectation value of the spin-antisymmetric 
orbital polarization P s = P-j- — P^ is non-zero, as well as 
total P = Pf + P_i- While P only breaks C4 symmetry, P s 
also break the time-reversal symmetry T (but preserves 
the product T x C\). It was proposed in the context of 
the nematic phase in Sr 3 Ru 2 07 that appearance of both 
P and P s could lead to an additional symmetry-allowed 
term in the Ginzburg-Landau free energy oc M(PP S ) = 
M(P 2 — P 2 ) [41]. Our analysis shows however that if 
it were the case, in the presence of orthorhombic lattice 
distortion P a oc 6 and hence M oc P 2 — P 2 oc <5 2 , instead 
of the experimentally observed M oc 5. 

In order to provide a more microscopic basis for the 
Ginzburg Landau theory Eq. (1), we studied the inter- 
play between orbital and spin physics within the frame- 
work of the Kugel-Khomskii model [42] . In the t 2g orbital 
basis, the Kugel-Khomskii model can be schematically 
written as follows [43, 44]: 

H = J 1 ]T [(S, • S l+X + l)T?Tf +A + (S { ■ S l+y + 1)t<Y^] 

(5) 



+ j 2 ]T (Si-s,-)- J d ^(s i -s j 



where Si = 1 is the magnetic moment localized on Fe 
d xz and d yz orbitals, and {r a ,T 6 } are the pseudospin 
operators that act in the orbital subspace of \xz) and 
\yz) states and depend on the directionality of the Fe-Fe 
bond. The last term describes the direct (ferromagnetic) 
exchange between Fe spins, proposed by R. Singh [18]. 

The key question is the occupation of the d xz and d yz 
orbitals. In Ref. [44], a similar model was derived for 
the configuration with total of 2 electrons in three i 2g 
orbitals, whose spins are aligned to form spin S = 1. 
However our GGA and GGA+Z7 calculation predict that 
the total number of electrons in the i 2g orbitals is about 
3.16 in the paramagnetic phase, with each orbital equally 
populated close to half filling. This agrees with the Fe 2+ 
multiplet structure calculated from LDA+DMFT [45]. 
Thus the Fe 2+ ground state appears to lie in the i 2g 
rather than t?, g sector, as was assumed in Ref. 44. Nev- 
ertheless, if one only considers \xz) and \yz) states, they 
can still be described by total spin S = 1 and Eq. (5) 
should hold. However the crucial difference from earlier 
treatments [18, 21, 44] is that n are not Ising variables, 
since both xz and yz orbitals are partially occupied. 

In what follows, we develop the mean-field treatment 
of the model (5), approximating the orbital average 



( T f T f+s> ~ (( T f ) 2 ) ~ n lz and similarly (t?7? +4 ) w n 2 yz . 
Eq. (5) can be decomposed into a sum of orbital and spin 
parts: H M f = E a + H a + P sp , where 

-Eo = Jr [< S * ' S *+* + !>»L + (Si ' S l+y + l)n 2 yz ] (6) 

i 

H Q = JiJ^Si • S l+£ + 1)t?7? +4 + (Si ■ S l+y + 1)tH 



b 

i+y 



Hsp — J]] [Jia Si • Si+x + Jib Si ■ Si +y ] + J 2 Sj-Sj(7) 

We see that the mean-field ground state energy takes on 
the form Eq <~ — 7(SiSj)n 2 1 , , justifying the choice of the 
biquadratic coupling in the Landau theory Eq. (1). 

It follows from H Q above that the strength of the in- 
tersite orbital coupling is proportional to (S, • + 1), 
i.e. the deviation of spins from the classical Neel state. 
This is a general feature of the Kugel-Khomskii model, 
well-known for the e g case [46, 47]. The classical Neel 
state with ordering wave- vector (ir, it) would have infinite 
orbital degeneracy corresponding to SU(2) pseudospin 
invariance. However for the case of the collinear (tt, 0) 
order observed in the pnictides, only the orbital correla- 
tions along the antiferromagnctic a-axis will be massless. 
These massless local orbital fluctuations will affect the 
spin sector, reducing the staggered moment [46, 47]. 

In the effective spin Hamiltonian (7), the exchange 
"constants" depend on the orbital polarization P: 



•ha = Ji(T?Ti +s .) - J d 

Jib = Jl(T~iT i+i ) - J d 



J 1 (n + P) 2 /A-J d (8) 
J!(n-P) 2 /4- J d , (9) 



where we denoted the total number of electrons in d xz 



and d yz orbitals by n 



. From our GGA and 



GGA+C7 calculations, we find n w 2.5 in the collinear or- 
dered state, in difference to Ref. 21 where it was assumed 
that orbital degrees of freedom were Ising and that n = 1. 

At the mean-field level, we can express (Si • Sj) 
through the ordered moment M, e.g.: (Sj-Sj+g) = ±M 2 , 
where the "+" sign corresponds to ordering wave-vector 
Q = (tt,0) and "-" sign to Q = (0,tt). The orbital po- 
larization is found by minimizing the ground state energy 
E in Eq. (6) with respect to n xz = (n + P)/2: 



;(i±m 2 ), p 



{- 



nM 2 , 
jiM 2 , 



for Q = (tt, 0) 
for Q = (0, tt) 



(10) 

The above two solutions are degenerate and the system 
chooses one of them by breaking the Z 2 Ising symmetry, 
as in the "spin nematic" scenario [17]. 

Within our ab initio GGA calculations, the ordered 
moment in the xz and yz orbitals M = 0.38. Using 
"gga = 2.54, Eq. (10) yields the mean-field value for 
orbital polarization Pmf — 0.37, which is significantly 
higher than Pqga = 0.19 from the ab initio calculation. 
This is not unexpected since orbital fluctuations will tend 
to lower the polarization from its mean field value. 
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We now return to the effective spin Hamiltonian, 
Eq. (7). At the mean field level, the anisotropy of the 
effective exchange "constants" in a6-plane becomes 

_ Jig ~ Jib = 2nP 

V ~ Jia + Jib ~ n 2 + P 2 - A\J d \/.h [ > 

Using the ARPES estimate [26] P w 0.18 (consistent 
with our Pgga = 0.19), n rts 2.5, and taking a real- 
istic estimate on |J<j| = 10 mcV, J\ = 50 meV [12], 
we arrive at the value r\ w 0.16. This would suggest 
that Jk, i=a 1.35Ji — J^, which is clearly antiferromag- 
netic, unless one assumes an unrealistically large value for 
Jd > J\. It has been proposed [18, 21] that orbital order 
coupled to Fe spins as in Eq. (7) could explain the sig- 
nificant anisotropy in the spin-wave dispersion observed 
by neutron scattering in CaFe2As2 and BaFe2As2 [5, 12]. 
Fitting neutron spectra to the linear spin- wave theory ac- 
tually results in a ferromagnetic Jib ~ —6 meV [5]. How- 
ever in Eq. (11), even for a large value of \ Jd\ = Ji/2, the 
upper bound on rj is <~ 0.21, and Ju, remains positive. 

Using the mean-field expression for orbital polariza- 
tion Eq. (10), we can express the parameters of the spin 
Hamiltonian in terms of the ordered moment as follows: 

Jla/6 = {^f J d ) ± + G(M 4 ) (12) 

What is striking is that the values of the exchange "con- 
stants" depend on the ordered moment. This has two 
consequences: first, the anisotropy in the spin-wave dis- 
persion grows in the ordered phase, in turn reinforc- 
ing the orbital polarization P cx M 2 . Secondly, the 
term ~ M 2 in Eq. (12) introduces quartic coupling be- 
tween neigbouring spins of the form (Si • Sj)(Sj • Sj), 
making the model qualitatively different from a simple 
anisotropic Heisenberg model. This strongly suggests 
that attempts to fit the neutron spin-wave spectra to 
a simple J\ a — Ju, — Ji model with constant exchange 
parameters [5, 12] are problematic, since they neglect to 
take quartic coupling between the spins into account. 

Present work suggests that orbital ordering, as seen 
in our ab initio calculations and in experiment [25, 26], 
is consistent with predictions of microscopic theory and 
explains a number of key experimental features, in partic- 
ular the linear correlation of ordered moment with lattice 
distortion, M ~ S. We predict that orbital polarization 
should scale as a square of magnetic moment close to 
Tjv, a feature that can be tested experimentally. We pro- 
pose that the spin- wave dispersion in the magnetically or- 
dered phase should be fitted with the moment-dependent, 
rather than constant, exchange parameters J\ a , Ju, given 
by Eq. (12). Although self-sufficient, orbital polarization 
can in principle co-exist with other proposed mechanisms 
for spontaneously broken £4 symmetry, such as spin-Ising 
ordering [14, 15] or Pomeranchuk instability [48]. 
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